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Abstract — A  diffusion-drift  treatment  of  electron  and  hole 
transport  in  macroscopic  graphene  is  presented.  The  various  ma¬ 
terial  response  functions  that  enter  the  theory  are  outlined  and, 
to  the  extent  possible,  specified  and  calibrated.  For  purposes  of 
illustration,  the  theory  is  applied  to  a  variety  of  situations  involving 
field-effect  devices  that  are  of  potential  technological  interest.  Both 
single  and  multilayer  graphene  are  discussed,  as  is  the  effect  of  the 
small  bandgaps  that  have  been  reported  for  graphene  on  SiC. 

Index  Terms — Diffusion  drift  (DD),  field-effect  transistors, 
graphene,  multilayer. 

I.  Introduction 

The  discovery  that  single-layer  graphene  is  a  robust 
material  that  can  be  isolated,  contacted,  and  tested  electri¬ 
cally  has  generated  much  excitement  [1],  [2].  Among  the  results 
obtained  to  date,  researchers  have  measured  mobilities  as  high 
as  250000  cm^/V  •  s  [3],  demonstrated  field-effect  transistor 
action  [1],  [4]-[6],  and  even  done  an  initial  characterization  of 
transport  in  graphene  nanoribbons  [7].  Most  recently,  intriguing 
evidence  has  been  obtained  for  the  existence  of  a  bandgap  of 
roughly  0.26  eV  for  single-layer  graphene  when  it  is  formed 
epitaxially  on  SiC  by  sublimation  of  silicon  at  high  temperature 
[8]-[10].  All  of  this  activity  is  clearly  promising,  but  of  course, 
whether  it  will  eventually  lead  to  a  practical  graphene  elec¬ 
tronics  technology  remains  to  be  seen.  Eor  this  paper,  we  take 
the  rapid  progress  that  has  been  made  as  sufficient  motivation 
for  framing  a  device-oriented  transport  theory  for  graphene 
that  may  one  day  serve  as  the  foundation  for  electronic  device 
design  and  optimization  tools. 

Theoretical  work  on  electron  transport  in  graphene  has 
generally  centered  on  understanding  ballistic  or  near-ballistic 
behaviors  and  on  the  exotic  electronic  properties  that  graphene 
can  exhibit  under  ideal  conditions,  including  when  restricted  in 
dimension  as  in  the  nanoribbons  [7].  In  this  paper,  we  focus 
instead  on  the  opposite  limit,  i.e.,  on  developing  a  device  mod¬ 
eling  approach  that  is  apropos  when  the  transport  is  dominated 
by  scattering.  This  interest  is  justified  in  part  by  the  fact  that 
except  for  the  technologically  uninteresting  case  of  suspended 
exfoliated  graphene  [3],  the  mobilities  measured  in  graphene 
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have  so  far  generally  been  far  from  ideal  [5],  [6].  It  seems  that, 
at  least  at  present,  the  all-surface  nature  of  graphene  makes 
it  particularly  susceptible  to  nearby  charges  and/or  localized 
“midgap”  states  [11].  While  this  issue  is  clearly  in  need  of 
attention,  it  does  imply  that  a  graphene  transport  description 
that  assumes  strong  scattering  will  be  relevant  at  least  for  most 
research  devices  of  today. 

The  conventional  transport  modeling  approach  with  strong 
scattering  is,  of  course,  the  diffusion-drift  (DD)  description 
[12].  This  theory,  which  was  developed  over  50  years  ago  and 
is  still  the  workhorse  for  practical  semiconductor  device  sim¬ 
ulation,  is  a  macroscopic  (or  continuum)  description  in  which 
the  electron  and  hole  populations  are  represented  by  continuous 
fluids  or  gases.  The  goal  of  this  paper  is  to  develop  a  similar 
DD  description  for  the  carrier  transport  in  graphene.^  Although 
ideal  graphene  is  a  semimetal  with  a  markedly  different  band 
structure  than  an  ordinary  semiconductor,  the  basics  of  its  DD 
description  are  not  much  different  from  the  conventional  theory. 
Eor  one  thing,  we  shall  regard  electrons  in  extended  states  above 
the  Dirac  point  as  forming  an  electron  continuum,  and  electrons 
absent  from  extended  states  below  the  Dirac  point  as  forming  a 
hole  continuum.  Moreover,  as  in  the  conventional  description, 
carrier  motion  in  the  graphene  will  be  assumed  to  be  determined 
by  the  net  effect  of  the  following  three  macroscopic  forces: 
1)  the  force  exerted  by  the  electrostatic  field  (drift);  2)  the 
statistical  pressure  force  arising  from  the  random  motions  of 
the  carriers  (diffusion);  and  3)  the  drag  force  generated  by 
scattering  that  retards  carrier  flow  through  the  lattice.  Within 
this  framework,  the  main  qualitative  difference  between  trans¬ 
port  in  graphene  and  in  conventional  semiconductors  is  the 
2-D  nature  of  single-layer  grapheme,  which  precludes  transport 
normal  to  the  layer.  (Eor  multilayer  graphene,  we  provide  here 
only  the  lowest  order  treatment  in  which  the  possibility  of 
transfer  of  carriers  between  layers  is  ignored,  and  the  model 
therefore  consists  merely  of  an  electrostatically  coupled  stack 
of  single  layers).  As  a  result,  the  DD  transport  equations  will 
be  partial  differential  equations  in  2-D.  In  all  other  respects, 
however,  the  DD  equations  for  graphene  are  familiar  with  the 
remaining  differences  arising  simply  from  graphene  being  a 
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different  “semiconductor”  with  a  different  band  structure  and 
transport  properties. 

Judging  from  its  conventional  counterpart,  the  DD  theory  of 
graphene  will  be  rigorously  justifiable  when  the  electron/hole 
mean  free  path  is  small  compared  to  the  device  size.  The  theory 
can  therefore  be  expected  to  apply  to  most  of  the  graphene 
test  devices  of  current  interest,  which  have  relatively  “large” 
geometries  and  are  composed  of  low-mobility  graphene  [5],  [6]. 
Beyond  this,  there  is  the  well-known  simplicity  and  robustness 
of  the  DD  approach  that  has  led  to  its  continued  utility  in 
electronics,  even  when  its  foundational  assumptions  no  longer 
seem  valid.  For  graphene,  it  seems  reasonable  to  suppose  that 
the  DD  description  will  remain  a  reliable  guide  even  for  higher 
mobility  (or  even  ballistic)  situations  since  it  will  still  incor¬ 
porate  the  correct  and  often-dominant  geometry,  band  structure 
effects,  and  electrostatics.  If  a  graphene  device  technology  is 
ever  realized,  the  DD  theory  of  graphene  may  well  prove  to 
be  a  useful  theory  and/or  phenomenology  for  practical  device 
design,  characterization,  and  optimization. 

II.  Equations  of  the  DD  Theory  of  Graphene 
A.  Balance  Laws 

Assuming  the  electrons  and  holes  in  graphene  can  be  re¬ 
garded  as  continuous  gases  with  well-defined  charge  and  mo¬ 
mentum  densities,  the  partial  differential  equations  imposed  on 
these  densities  by  charge  conservation  will  take  the  following 
form: 

V  •  Jn  =  —  where  =  —qnwn 

V  •  Jp  =  —  i?,  where  =  qpw^.  (1) 

While  looking  like  the  ordinary  3-D  equations,  because  the 
carriers  are  restricted  to  be  within  the  graphene  layer  (whether 
it  is  planar  or  not),  these  are  2-D  equations  where  n  and  p  are 
the  areal  densities  (in  units  per  square  centimeter),  and  Vp 
are  the  gas  velocity  vectors  (in  centimeters  per  second)  lying  in 
the  2-D  surface,  and  Jp  are  the  current  density  vectors  (in 
amperes  per  centimeter),  R  is  the  recombination/generation  rate 
(in  coulombs  per  square  centimeter  second),  and  the  differential 
operators  are  surface  divergences. 

Scattering-dominated  flow  means  that  electron  and  hole  iner¬ 
tia  can  be  neglected  and  the  momentum  balance  equations  for 
the  electron  and  hole  gases  then  reduce  to  force  balances  among 
the  three  basic  forces  mentioned  in  Section  I.  Formulating  the 
pressure  force  within  each  carrier  gas  in  terms  of  the  chemical 
potentials  pn  and  Pp,  and  assuming  the  drag  forces  per  charge 
are  proportional  to  the  gas  velocities  (with  the  usual  mobilities 
Pn  and  Pp  as  proportionality  constants),  we  can  then  recast  the 
force  balances  as  the  following  DD  equations: 

J„  =  gn/i„V(v2„  -  tA)  3p  = -qnidpV{(fip  +  ip)  (2) 

where  is  the  electrostatic  potential.  As  with  (1),  the  equations 
in  (2)  are  familiar  in  form  but  are,  again,  2-D  equations  with 
the  derivatives  being  surface  gradients  within  the  graphene 
layer.  When  more  than  a  single  layer  of  graphene  is  present. 


each  extra  layer  would  be  described  by  additional  equations 
analogous  to  those  in  (1)  and  (2),  and  with  the  rate  R  in  (1) 
possibly  also  incorporating  interlayer  transfer/recombination 
processes. 

Finally,  because  the  graphene  layers  are  generally  embedded 
within  a  3-D  device  structure  that  includes  insulators  and  metal 
contacts,  the  2-D  equations  describing  transport  in  the  graphene 
must  be  coupled  to  a  3-D  electrostatic  equation,  i.e., 

V  •  (edV?/))  =  -gA^fixed  (3) 

where  A/fixed  and  Sd  are  the  position-dependent  fixed  vol¬ 
umetric  charge  density  and  permittivity,  respectively,  in  the 
surrounding  materials. 

In  addition  to  the  differential  equations,  the  basic  theory 
must  also  include  a  consistent  set  of  boundary  conditions.  In 
general,  these  are  expressions  of  the  fundamental  balance  laws 
already  discussed  and  are  derived  by  taking  appropriate  limits 
of  integral  versions  of  them.  However,  because  the  derivations 
and  resulting  conditions  are  mostly  familiar  [12],  we  omit 
further  discussion,  except  to  note  that  the  crucial  electrostatic 
condition  derived  from  Gauss’  law  that  relates  the  electric 
displacements  across  the  (or  each)  graphene  layer  is 

n  ■  {£d2^2  -  Sdi'Ei)  =  q{p  -  n  +  N)  (4) 

where  the  subscripts  indicate  values  on  either  side  of  the 
graphene  layer,  n  is  the  normal  vector  pointing  from  the  “1“  to 
the  “2”  side  of  the  layer,  and  N  is  the  effective  surface  charge 
due  to  any  dopants,  whether  substitutional  or  adsorbed  [13],  or 
fixed  charges  at  the  layer/interface. 

B.  Material  Response  Functions 

The  equations  in  Section  II-A  are  expressions  of  the  general 
conservation  laws  of  classical  physics  as  applied  to  the  electron 
and  hole  continua  in  graphene,  assuming  strong  scattering 
(i.e.,  scattering  forces  dominating  over  inertial  forces).  In  large 
part,  they  represent  the  physical  content  of  the  DD  theory. 
However,  as  is  usual  in  classical  field  theory,  these  equations 
are  not  complete  in  both  physical  and  mathematical  senses. 
To  complete  the  DD  theory  of  graphene,  one  must  particu¬ 
larize  various  material  response  functions  that  introduce  into 
the  general  equations  the  specific  properties  of  graphene.  It 
is  important  to  note  that  because  the  growth  techniques  and 
properties  of  “device-grade”  graphene  are  yet  to  be  defined 
(and,  indeed,  may  never  be),  optimal  expressions  and  coeffi¬ 
cient  values  for  these  material  response  functions  cannot  be 
established  yet.  Therefore,  the  functions  presented  herein,  and 
their  consequences  discussed  in  subsequent  sections,  should 
be  regarded  merely  as  illustrative.  As  knowledge  of  graphene 
transport  advances,  effort  should  focus  on  better  defining  these 
response  functions  and  on  measuring  their  coefficient  values. 
Of  course,  in  the  process,  the  range  and  limitations  of  the  DD 
framework  itself  for  real  graphene  will  likely  also  be  revealed. 

The  functional  dependences  of  the  material  response  func¬ 
tions  that  appear  in  the  DD  equations  are  largely  set  by  thermo¬ 
dynamic  considerations.  This  connection  will  not  be  discussed 
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here  [14];  instead,  we  note  only  the  following  primary  func¬ 
tional  dependences: 


D  =D(E,T) 

(5a) 

9>n  =‘Pn{n,T) 

‘Pp  =  ‘Pp{P,T) 

(5b) 

II 

II 

(5c) 

R  =R{n,p,  T). 

(5d) 

Of  these  material  response  functions  (also  known  as  consti¬ 
tutive  equations),  the  ones  for  which  specific  expressions  are 
most  easily  developed  are  those  in  (5  a)  and  (5b)  since  they  are 
nontrivial  under  equilibrium  conditions  and  do  not  involve  the 
nonequilibrium  quantities  and  Vp. 

Equation  (5a)  is,  of  course,  already  familiar  from  classical 
electrostatics,  and  its  usual  linear  form  D  =  has  already 
been  incorporated  in  (3)  and  (4).  With  respect  to  the  “equations 
of  state”  in  (5b),  as  in  ordinary  DD  theory,  their  equilibrium 
nature  allows  us  to  utilize  the  powerful  tools  of  statistical 
mechanics.  In  particular,  we  start  with  the  general  formula  for 
the  density  of  conduction  electrons  (with  the  analysis  for  holes 
being  completely  analogous)  in  uniform  equilibria,  i.e., 

oo 

—  f _ 9c{E)dE _ 

^  ^  “  y  1  +  exp  [(E  -  EF)kBT] 

0 

where  gc{E)  is  the  2-D  density  of  states  in  the  conduction  band. 
In  ideal  single-layer  graphene,  the  expression  for  gc{E)  is  [1], 
[15] 

2E 

where  the  energy  is  measured  from  the  Dirac  point  and  the 
Dirac  electron  velocity  c  in  graphene  is  approximately  10^  cm/s 
[1].^  When  graphene  has  a  bandgap  Eq  (again,  as  a  result  of 
its  association  with  SiC  [8]-[10]),  we  employ  the  following 
approximation  to  the  conduction  band  structure: 

E{k)  =  +rflk^c^ 


where  the  energy  is  now  measured  from  midgap  and  it  is 
assumed  that  both  valleys  remain  degenerate  as  suggested  by 
[9] .  It  is  then  readily  shown  that  the  conduction  band  density  of 
states  is  given  by 


9c{E)  = 


2 

7rJi?c^ 


E  <  Eq/2 
otherwise 


(7b) 


which  reduces  to  (7a)  when  Eg  vanishes. 


^In  multilayer  graphene,  interlayer  coupling  can  modify  the  band  struc¬ 
ture,  as  has  been  observed  in  both  photoemission  [T.  Ohta,  A.  Bostwick, 
J.L.  McChesney,  T.  Seyller,  K.  Horn,  and  E.  Rotenberg,  “Interlayer  interac¬ 
tion  and  electronic  screening  in  multilayer  graphene  investigated  with  angle- 
resolved  photoemission  spectroscopy,”  Phys.  Rev.  Lett.  98,  206802  (2007)]  and 
cyclotron  resonance  [E.A.  Henrikson,  Z.  Jiang,  L.-C.  Tung,  M.E.  Schwartz, 
M.  Takita,  Y.-J.  Wang,  R  Kim,  and  H.L.  Stormer,  “Cyclotron  resonance  in 
bilayer  graphene,”  Phys.  Rev.  Lett.  100,  087403  (2008)].  This  effect  is  ignored 
in  this  paper  both  for  simplicity  and  because  its  impact  is  primarily  at  energies 
greater  than  the  coupling  energy  estimated  to  be  0.39  eV. 


10^  10®  10®  10’®  10^^  10’^  10’® 
N[j(cm-®) 


0.7 
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0,5 


0.4  m 

0.3  ^ 
0.2  < 


0.1 

0 


-0.1 


Fig.  I.  Electron  and  hole  densities  and  the  Fermi  energy  versus  the  doping 
density  in  single-layer  graphene  both  with  the  ideal  zero-gap  band  structure 
and  with  a  bandgap  of  0.26  eV. 


Fig.  2.  Intrinsic  density  and  the  coeffcient  Nq  in  (8)  versus  bandgap  in  single¬ 
layer  graphene. 


By  varying  the  carrier  density  in  the  graphene  via  doping 
{N  =  Nd),  we  can  study  the  relationship  between  density  and 
Fermi  level  implied  by  the  foregoing  equations.  Inserting  (7b) 
into  (6),  with  Eq  =  0  eV  or  0.26  eV,  we  obtain  numerically 
the  results  shown  in  Fig.  1.  Of  course,  having  a  bandgap  is 
desirable  from  a  device  perspective,  and  this  may  be  seen  in  the 
reduction  in  the  intrinsic  (neutral)  density,  as  plotted  in  Fig.  1. 
More  explicitly,  in  Fig.  2,  we  show  the  direct  dependence  of  the 
intrinsic  density  on  bandgap. 

For  use  in  the  DD  equations,  we  need  an  expression  for  (5b) i 
(and  analogously  for  (5b)2),  and  the  following  approximate 
form  is  easily  reached  by  numerical  inversion: 


T)  =  ksT 


In 


n/Ngr  \ 
\EEG/2kBT)^ 


n/Ngr  Y1 
lEN^/n)  I 


(8) 
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Fig.  3.  Variation  of  the  electron  density  as  a  function  of  either  the  exact 
chemical  potential  (lines)  or  the  approximation  of  (8)  for  ideal  graphene 
(diamonds)  and  graphene  with  a  bandgap  of  0.26  eV  (squares).  The  relative 
errors  in  each  case  are  also  shown. 


for  the  material  response  functions  in  (5c)  and  (5d)  because 
of  their  nonequilibrium  nature,  and  one  must  instead  rely 
on  complicated  experimental  measurements  and/or  theoretical 
calculations.  Unfortunately,  neither  experiments  (not  being  suf¬ 
ficiently  comprehensive  or  repeatable)  nor  theory  (not  being 
sufficiently  realistic)  are  presently  adequate  for  this  purpose. 
Therefore,  in  this  paper,  we  consider  only  simple  choices  for 
(5c)  and  (5d),  with  the  goal  of  illustrating  the  theory  and  with 
no  claim  that  these  choices  are  actually  representative  of  real 
graphene.  Presumably,  as  research  on  graphene  progresses, 
either  these  forms  will  be  confirmed  or  better  choices  can  be 
developed. 

With  respect  to  the  mobilities,  two  specific  choices  are 
considered.  In  one,  we  simply  assume  constant  mobilities,  for 
definiteness  taking  them  to  be  /i^o  =  MpO  =  500  cm^/V  •  s. 
This  selection  is  particularly  useful  when  the  impact  of  other 
constitutive  choices  is  being  investigated.  The  other  mobility 
model  follows  the  conventional  approach  of  substituting  elec¬ 
tric  field  for  velocity  and  incorporates  by  an  inverse  sum  a 
dependence  on  the  component  of  the  electric  field  normal  to 
the  graphene  layer  {E_i),  i.e.. 


Fig.  4.  Enhancement  factor  for  diffusion  in  graphene  due  to  its  bandstructure. 

where  r  =  0.58,  R  =  0.69,  Nq  varies  with  Eq  as  shown 
in  Fig.  2,  and  the  effective  density  of  states  is  Ngr  = 
2{kBT/hcY /tt  which  at  room  temperature  is  approximately 
equal  to  10^^  cm“^.  As  seen  in  Fig.  3,  (8)  accurately  de¬ 
scribes  the  electron  gas  for  Eg  =  0  and  0.26  eV  for  n  < 
5  X  10^^  cm“^.  Due  to  the  symmetry  of  the  graphene  band 
structure,  one  can  immediately  write  an  equation  of  state  for 
holes  that  is  directly  analogous  to  (8).  Finally,  from  (8),  we 
can  derive  an  expression  for  the  diffusion  “constant”  using 
the  formula  =  find^n/dn.  A  plot  of  the  diffusion  “en¬ 
hancement”  factor  Dn/iin  is  shown  in  Fig.  4,  and  as  one 
can  see,  the  band  structure  of  graphene  greatly  amplifies  the 
diffusivity  of  electrons  as  the  density  increases.  We  note  that 
the  widening  of  screening  layers  that  would  be  associated  with 
this  enhanced  diffusion  may  well  have  consequences  for  the 
ultimate  scalability  of  graphene  devices. 

As  in  conventional  DD  theory,  it  is  not  possible  to  be 
nearly  as  definitive  when  it  comes  to  developing  expressions 
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where  E^,  Ep,  7^,  and  7^  are  constants.  In  Si  MOSFETs,  such 
electric  field  dependences  are  commonly  included  to  account 
for  interface  roughness  scattering.  For  graphene,  however, 
should  such  dependences  actually  prove  important  for  describ¬ 
ing  the  transport,  they  seem  more  likely  to  arise  from  scattering 
by  fixed  charges.  This  view  is  supported  by  the  simulations  of 
[16],  [17],  which  suggest  that  small  changes  in  the  distance 
between  the  channel  charge  and  the  fixed  charge  (on  the  order 
of  1-2  A)  could  significantly  impact  the  scattering. 

Because  of  graphene’s  zero/small  bandgap,  the 
generation-recombination  term  R  is  undoubtedly  almost 
always  important  in  graphene  devices  when  they  are  driven 
out  of  equilibrium.  Many  mechanisms  may  be  involved,  but 
favoring  simplicity,  we  consider  just  two — one  associated  with 
thermal  generation  and  the  other  with  band-to-band  tunneling. 
The  specific  expressions  we  assume  are 
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where  t,  a,  (3,  and  7  are  constants,  and  the  two  terms  combine 
SLS  R  =  Rth  +  ^tun-  Equation  (10a)  is  arbitrarily  assumed  to 
be  of  the  Shockley-Read-Hall  form  [12],  and  the  band-to-band 
tunneling  expression  in  (10b)  is  also  based  on  conventional 
expressions  [18].  The  latter’s  dependence  on  the  difference 
in  quasi-Eermi  levels  ensures  that  it  vanishes  in  equilibrium 
and  that  it  has  the  proper  dependence  on  occupancy  (e.g.,  net 
recombination  when  Ep  >  Ep).  Based  on  the  conventional 
models,  we  assert  dependences  on  Eg  as  follows: 
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Fig.  5.  Depictions  of  graphene  field  effect  transistors  simulated  in  this  paper, 
(a)  Model  for  the  back-gated  device  studied  in  [5]  and  [6].  (b)  Model  for  a 
front-gated  device  such  as  those  formed  on  SiC  substrates. 
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where  ctp  and  f3o  are  constants,  and  A  is  a  small  constant  energy 
ensuring  that  a  does  not  become  unbounded  as  Eg  0. 

Finally,  we  observe  that  in  the  treatment  of  multilayer 
graphene,  interlayer  electron  transfers  would  enter  the  the¬ 
ory  through  an  additional  generation-recombination  term 
^interlayer  (iD  addition  to  the  possible  band  structure  effects 
mentioned  in  [16]).  Such  transfers  are  ignored  in  this  paper  both 
for  reasons  of  simplicity  and  because  we  believe  their  role  to  be 
negligible  for  the  situations  considered  in  this  paper  in  which 
the  highly  conductive  layers  are  relatively  uniform  and  always 
share  common  contacts. 

III.  DD  Simulations  of  Graphene 
Device  Phenomena 

A.  Field  Effect 

To  examine  the  DD  treatment  of  the  field  effect  in  graphene, 
we  model  the  particular  FET  studied  in  [5]  and  [6],  shown 
schematically  in  Fig.  5(a).  In  this  device,  the  Si  substrate  is 
doped  and  contacted  so  that  it  can  function  as  a  back  gate.  The 
experimentally  measured  drain  current  at  low  drain  voltage  is 
plotted  versus  the  applied  gate  field  in  Fig.  6  (squares)  [5],  [6]. 
Clearly,  a  significant  field  effect  is  observed,  but  not  surpris¬ 
ingly  given  the  lack  of  a  bandgap,  the  device  displays  poor  turn¬ 
off  characteristics.  In  fact,  the  lack  of  rectifying  contacts  results 
in  bipolar  conduction  with  the  current  merely  going  through 
a  minimum  as  the  channel  switches  from  electron-  to  hole- 
rich,  as  shown  in  the  inset  in  Fig.  6.  Also  noteworthy  are  the 
asymmetries  in  this  characteristic  both  in  voltage  and  in  current. 

In  modeling  the  device  in  Fig.  5(a),  for  simplicity,  we  ignore 
possible  nonuniformities  in  the  width  direction  (for  now)  and 
apply  the  DD  equations  in  reduced  dimension,  specifically 
solving  (1)  and  (2)  in  1-D  and  (3)  in  2-D.  These  equations  can 
be  solved  only  numerically,  and  for  this  purpose,  we  transform 
to  Slotboom  variables  to  improve  conditioning  [12]  and  then 
utilize  the  finite-element  method  on  a  nonuniform  unstruc¬ 
tured  grid.  In  general,  the  simulations  always  produce  weak 
field-effect  characteristics  like  those  in  Fig.  6,  showing  the 
characteristic  bipolar  conduction  and  a  high  OFF-state  current. 
However,  the  quantitative  details  including  the  voltage/current 
asymmetries  depend  on  the  specific  choices  for  the  dopants. 


0 . 111. . . . .  I  ii 

-3  -2-10123 

Electric  Field  (MV/cm) 

Fig.  6.  Experimental  data  of  [5],  [6]  and  fit  using  a  I-D/2-D  DD  simulation,  as 
discussed  in  the  text.  The  inset  shows  the  simulated  electron  and  hole  currents 
in  the  device. 

the  fixed  charge,  the  contact  properties,  the  mobility  models, 
and  the  generation-recombination  models.  Inadequate  knowl¬ 
edge  of  the  device  prevents  us  from  being  definitive  in  these 
choices  and  from  identifying  which  of  them  is  most  important. 
Therefore,  with  illustration  in  mind,  we  simply  assume  that  the 
observed  threshold  voltage  shift  is  due  to  fixed  charge  at  the 
graphene-Si02  interface  (with  N  =  —10^^  cm“^  ^  AVr  = 
0.4)  and,  following  [5]  and  [6],  that  the  mobilities  vary  with 
the  vertical  electric  field  as  given  in  (9).  Furthermore,  in  order 
to  provide  a  well-defined  contact,  we  assume  that  the  graphene 
in  the  vicinity  of  the  contacts  is  doped  n-type  {N  =  Nd)-  That 
such  assumptions  allow  DD  theory  to  fit  the  data  is  shown  in 
Fig.  6  (line),  where  we  have  assumed  //nO  =  1690  cm^/V  •  s, 
7n  =  0.5,  En  =  0.5  MV/cm,  //po  =  3360  cm^/V  •  s,  7n  = 
0.25,  and  Ep  =  0.5  MV/cm.  Again,  it  is  not  known  whether 
these  strong  electric  field  dependences  are  actually  responsible 
for  the  observed  behavior,  but  it  is  certainly  plausible. 

B.  Effect  of  Bandgap 

As  noted  in  Section  I,  a  bandgap  is  believed  to  be  generated 
in  graphene  when  it  is  formed  epitaxially  via  high-temperature 
Si  sublimation  of  SiC  [8]-[10].  Having  a  bandgap  is,  of  course, 
attractive  because  it  is  its  lack  that  produces  the  small  on/off 
current  ratio  and  bipolar  conduction  seen  in  Fig.  6.  That 
graphene  on  SiC  is  believed  to  have  a  bandgap  as  large  as 
0.26  eV  [8]  makes  this  prospect  particularly  intriguing,  given 
that  this  is  roughly  the  same  as  the  bandgap  of  InSb,  which 
is  a  material  currently  considered  to  be  a  legitimate  candidate 
for  future  low-power  digital  electronics  [19].  To  illustrate  the 
bandgap  effect,  we  again  study  transport  in  an  FET  structure 
with  the  geometry  now  being  that  of  Fig.  5(b),  where  the 
substrate  is  SiC  and  the  frontside  has  an  insulated  gate  situated 
atop  the  graphene.  Obviously,  a  major  effect  of  the  bandgap  is 
on  generation-recombination,  and  we  include  this  effect  in  the 
simulations  using  (10a)  and  (10b)  with  (11).  Finally,  in  order 
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Gate  Voltage  (V) 

Fig.  7.  DD  simulation  of  the  effect  of  a  bandgap  on  the  I-V  characteristics 
of  a  single-layer  graphene  FET,  with  L  =  0.5  /^m  and  Vd  =50  mV. 


Position  (pm) 

Fig.  8.  Current  and  recombination  rate  as  a  function  of  source-drain  bias  in 
the  FET  simulated  in  Fig.  7  (no  bandgap)  and  biased  so  that  the  channel  is 
hole-rich  (Vq  =  —1  V).  Current  flows  are  due  to  generation-recombination 
processes,  particularly  band-to-band  tunneling. 

not  to  obscure  the  bandgap  effect  with  extraneous  physics,  we 
assume  constant  mobilities  in  these  simulations. 

The  simulated  transfer  characteristics  of  the  graphene  FET 
are  plotted  in  Fig.  7  for  an  applied  drain  voltage  of  50  mV 
and  with  the  bandgap  treated  as  a  parameter.  Under  a  positive 
gate  bias,  the  channel  region  is  electron  rich,  current  flows 
easily,  and  the  bandgap  has  little  effect.  As  the  gate  voltage  is 
reduced,  the  current  drops  due  to  the  field  effect.  However,  as 
the  bias  goes  negative,  the  effect  of  the  bandgap  becomes  more 
pronounced.  If  the  bandgap  is  nonexistent  or  small,  then  as  seen 
earlier,  strong  band-to-band  tunneling  provides  good  electrical 
contact  to  the  hole-rich  channel  and  allows  a  bipolar  current  to 
flow.  This  action  is  depicted  in  profile  in  Fig.  8  for  the  case  with 
zero  bandgap.  However,  when  the  bandgap  becomes  larger. 


Drain  Voltage  (V) 

Fig.  9.  Drain  characteristics  with  and  without  carrier  generation- 
recombination  and  showing  the  effect  of  the  bandgap  in  reducing 
multiplication. 

the  decline  in  band-to-band  tunneling  strongly  suppresses  the 
current  under  negative  bias,  as  seen  in  Fig.  7.  Indeed,  for 
bandgaps  greater  than  about  0.26  eV,  the  FET  exhibits  good 
turn-off  characteristics  with  increasingly  large  on/off  current 
ratios.  In  these  cases,  the  off-current  seen  in  Fig.  7  is  largely  due 
to  the  abruptness  of  the  source-drain  contact  doping  profile. 

One  additional  effect  of  the  small  bandgap  is  carrier  multi¬ 
plication  under  a  high  drain  bias.  For  purposes  of  illustration, 
we  assume  that  this  phenomenon  occurs  solely  via  the  genera¬ 
tion  processes  included  in  (10a)  and,  more  particularly,  (10b). 
Simulation  results  comparing  the  currents  with  and  without 
generation  included  and  for  several  bandgaps  are  shown  in 
Fig.  9.  Not  surprisingly,  having  a  bandgap  strongly  suppresses 
the  multiplication  effect  as  well,  particularly  at  a  low  drain 
voltage. 

C.  Effect  of  Multilayers 

To  illustrate  the  DD  simulation  of  transport  in  multilayer 
graphene,  we  again  assume  an  FET  structure  like  that  of 
Eig.  5(b),  but  now  with  the  current  flowing  in  a  stack  of  four 
graphene  layers,  all  identically  contacted  at  the  source  and 
drain  and  all  with  identical  constant  mobilities.  As  noted  earlier, 
we  assume  in  this  paper  that  the  transfer  of  carriers  between 
layers  is  impossible,  and  to  model  the  electrostatic  interactions 
explicitly,  we  assume  that  the  layers  are  separated  by  a  very 
small  distance,  with  the  condition  in  (4)  applied  across  each 
individual  layer. ^  The  separation  distance  was  assumed  to  be 
either  0.35  nm  (the  actual  value)  or  1.0  nm  (a  relaxed  value 
that  eases  mesh  requirements),  and  not  surprisingly,  both  values 
give  essentially  the  same  results.  Simulations  are  carried  out 
just  as  before,  and  in  Eig.  10,  we  exhibit  the  resulting  gate 
characteristics,  plotting  the  current  carried  in  each  individual 
layer  for  the  case  when  all  layers  have  zero  bandgap.  The 

^This  calculation  is  therefore  the  macroscopic  analog  of  that  of  P.B.  Visscher 
and  L.M.  Falicov,  “Dielectric  screening  in  a  layered  electron  gas,”  Phys.  Rev.  B 
3,  2541  (1971). 
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Fig.  10.  DD  simulation  of  the  current  contributions  from  the  individual  layers 
in  a  four-layer  graphene  FET  as  a  function  of  gate  voltage.  In  this  simulation, 
all  the  layers  are  assumed  to  have  zero  bandgap. 

behavior  is  readily  understood.  For  gate  voltages  near  threshold 
(in  this  case,  Vt  ^  —0.4  V),  the  intrinsic  density  in  each  layer 
leads  each  to  contribute  roughly  equally,  while  at  a  higher 
bias  (either  positive  or  negative),  the  accumulation  of  electrons 
or  holes  in  the  layer  closest  to  the  gate  causes  it  to  screen 
out  the  other  layers  and  to  be  increasingly  dominant.  From 
the  perspective  of  materials  characterization,  this  observation 
may  be  of  value  in  that  it  means  that  the  single  (topmost) 
layer  conduction  can  be  monitored,  even  when  multilayers  are 
present  and  simultaneously  contacted. 

As  noted  earlier,  studies  of  graphene  formed  on  SiC  suggest 
that  an  energy  bandgap  is  opened  up  due  to  the  presence  of 
the  SiC  [8].  These  studies  further  indicate  that  for  multilayer 
graphene  the  induced  bandgap  decreases  in  the  layers  more  dis¬ 
tant  from  the  SiC  with  the  approximate  values  for  the  bandgaps 
reported  in  [8]  as  follows:  first  layer,  0.26  eV;  second  layer, 
0.14  eV;  third  layer:  0.066  eV;  and  fourth  layer:  0.02  eV. 
As  a  second  multilayer  calculation,  we  compare  in  Fig.  11 
the  transfer  characteristics  for  four-layer  graphene  with  no 
bandgap  (labeled  “ideal”)  with  two  additional  simulated  curves 
(labeled  “normal”  and  “inverted”)  in  which  the  realistic  variable 
bandgaps  on  SiC  have  been  included.  The  “normal”  curve  refers 
to  a  device  like  that  of  Fig.  5(b)  in  which  the  SiC  is  the 
device  substrate  and  the  gate  barrier  is  a  deposited  material 
such  as  Si02,  whereas  the  “inverted”  case  has  the  SiC  on 
which  the  graphene  has  formed  as  the  barrier  layer  and  the 
device  substrate  being  a  deposited  or  wafer-bonded  material.  To 
our  knowledge,  the  latter  “inverted”  structure  has  not  yet  been 
explored  experimentally,  but  as  Fig.  1 1  shows,  the  situation  is 
interesting  because  it  shows  a  potential  for  significant  turn-off 
despite  the  graphene  stack  including  zero-bandgap  layers.  To 
understand  how  this  can  be,  we  observe  first  that  under  a  posi¬ 
tive  bias,  the  currents  for  all  cases  are  nearly  the  same  because, 
as  seen  earlier,  the  bandgap  has  little  effect  when  the  conduc¬ 
tion  is  entirely  by  majority  carriers  (in  this  case,  electrons). 
However,  as  we  saw  in  Figs.  7  and  8,  under  a  negative  bias. 


Gate  Voltage  (V) 


Fig.  II.  DD  simulation  of  the  total  current  in  a  four-layer  zero-gap  graphene 
FET  as  a  function  of  gate  voltage. 

the  current  flow  depends  on  strong  generation-recombination, 
which  is  suppressed  when  a  significant  bandgap  is  present.  In 
the  “normal”  case,  the  layer  with  the  biggest  bandgap  is,  farthest 
from  the  gate,  and,  therefore  the  current  suppression  is  of  least 
consequence  since  this  carries  minimal  current  anyway  due  to 
the  screening  (as  shown  in  Fig.  10).  On  the  other  hand,  in  the 
“inverted”  case,  the  bandgap  is  in  the  layer  closest  to  the  gate, 
and  therefore,  suppression  of  its  current  will  have  a  big  effect 
since  this  layer  is  the  dominant  carrier  of  current  (except  near 
threshold)  because  of  screening. 

Another  way  of  looking  at  the  results  in  Fig.  1 1  is  to  compute 
the  ratios  of  the  currents  in  the  ON  and  OFF  states.  If  we  assume 
that  Vdd  =  1-5  V  and  use  a  “2:1”  rule  [21]  to  define  Vqn  = 
Vt  +  1  V  and  Vqff  =  Vt  —  0.5  V,  then  for  Vd  =  50  mV,  we 
obtain  the  following  /qn/^off  ratios:  zero  gap,  3.9;  “normal,” 
3;  and  “inverted,”  13.  Obviously,  from  a  device  perspective,  the 
“inverted”  case  is  best,  and  while  not  ideal,  it  might  provide 
the  simplest  route  to  “useful”  graphene  devices.  In  addition, 
as  noted  before,  it  could  also  serve  as  a  means  for  better 
characterizing  graphene-on- SiC  material  in  nonideal  multilayer 
samples. 

D.  Effect  of  Lateral  Nonuniformities 

In  ultrasmall-geometry  FETs,  an  important  technological 
issue  is  the  impact  of  random  dopants/charges  in  producing  sig¬ 
nificant  device-to-device  variations  in  threshold  voltage.  Given 
the  all-surface  nature  of  graphene,  one  might  expect  such  an 
issue  to  be  even  more  important  in  graphene  devices  [16], 
[17].  In  this  paper,  we  do  not  explore  this  issue  at  all;  instead, 
we  simply  use  it  as  a  means  of  illustrating  the  DD  theory  of 
graphene  in  a  3-D  situation.  To  this  end,  we  insert  a  random 
distribution  of  fixed  charges  into  a  3-D  version  of  the  device 
shown  in  Fig.  5(a)  and  then  simulate  it  by  solving  the  equations 
of  our  DD  theory,  namely,  (1)  and  (2)  solved  in  2-D  and  (3) 
solved  in  3-D.  A  result  from  these  simulations  is  plotted  in 
Fig.  12,  where  we  show  the  electron  density  as  a  function  of 
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Fig.  12.  Electron  density  in  a  channel  with  randomly  positioned  fixed  charges 
in  a  2-D/3-D  DD  graphene  simulation.  The  region  shown  is  50  nm  x  50  nm  in 
size. 

position  in  the  channel  of  a  small-geometry  FET.  Perhaps  most 
noteworthy  is  the  small  size  of  the  perturbations  induced  by 
the  charges  as  a  result  of  their  being  in  very  close  proximity 
to  the  2-D  graphene  layer.  This  suggests  that  the  effect  of 
lateral  nonuniformities  in  graphene  might  actually  be  less  than 
in  conventional  FETs.  It  is  also  reminiscent  of  the  experimental 
observations  of  Martin  [21]. 

IV.  Final  Remarks 

This  paper  has  introduced  a  DD  theory  that  describes  elec¬ 
tron  and  hole  transport  in  single  and  multilayer  graphene.  The 
main  difference  in  formulation  from  ordinary  DD  theory  stems 
from  the  2-D  nature  of  macroscopic  grapheme,  which  causes 
the  new  theory  to  be  a  hybrid  of  2-D  transport  equations  and 
3-D  electrostatics.  The  other  differences  from  conventional 
theory  are  more  quantitative  in  nature  and  largely  result  from 
the  different  band  structure  of  graphene.  In  other  words,  in  this 
respect,  graphene  is  merely  a  different  “semiconductor.”  The 
graphene  properties  enter  the  theory  through  material  response 
functions  that  we  described  first  in  a  general  fashion  and  then  by 
proposing  some  specific  forms.  Given  the  presently  incomplete 
state  of  knowledge  about  graphene,  the  best  understood  of  these 
functions  are  the  equations  of  state  for  electrons  and  holes,  and 
we  provide  expressions  both  for  the  case  of  ideal  graphene 
with  zero  bandgap  and  for  graphene  that  has  a  bandgap  as  is 
believed  to  occur  when  it  is  formed  on  SiC.  New  and  more 
tentative  constitutive  models  for  the  mobilities  and  for  gen¬ 
eration/recombination  by  thermal  means  and  by  band-to-band 
tunneling  were  also  proposed.  The  DD  theory  of  graphene  was 
then  illustrated  using  numerical  simulations  that  exhibited  its 
treatments  of  field-effect  situations,  the  effect  of  bandgap,  the 
effect  of  multilayers,  and  the  effect  of  lateral  nonuniformities. 

Altogether,  the  DD  theory  of  graphene  is  quite  familiar,  and 
given  the  close  analogy  to  ordinary  DD  theory,  it  is  obviously 
well  suited  to  device  modeling,  particularly  for  fast-turnaround 
engineering  applications.  In  this  use,  it  is  important  that  one 
be  aware  of  the  chief  foundational  assumption  of  the  theory, 
which  is  that  the  carrier  mean  free  path  is  presumed  to  be 
“short”  compared  to  the  device  size.  This  means  that  it  is 
best  applied  to  larger  device  geometries  and  to  relatively  low 


mobility  graphene  like  most  present-day  material.  Beyond  this, 
however,  experience  with  conventional  DD  theory  suggests 
that  the  graphene  theory  will  be  robust  and  likely  useful  well 
beyond  its  supposed  range  of  applicability,  particularly  when 
used  as  an  engineering  phenomenology.  Additionally,  just  as 
conventional  DD  theory  is  often  extended  to  apply  to  transport 
in  confined  geometries  (e.g.,  quantum  wells),  so  might  our 
graphene  transport  theory  be  generalized  to  treat  transport  in 
the  quantum-confined  regime,  e.g.,  to  graphene  nanoribbons. 
Of  course,  whether  any  of  this  is  ever  important  depends  on  the 
extent  to  which  graphene  matures  as  an  electronic  material  in 
the  coming  years. 
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